Point defects and clustering in uranium dioxide by LSDA+U calculations 
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A comprehensive investigation on point defects and their clustering behavior in nonstoichiometric 
uranium dioxide U02±a: is carried out using LSDA+U method based on density functional theory. 
Accurate energetic information and charge transfers available so far are obtained. With these en- 
ergies that have improved more than 50% over that of pure GGA and LDA, we show the density 
functional theory predicts the predominance of oxygen defects over uranium ones at any compo- 
sitions, which is possible only after treated the localized 5/ electrons properly. Calculations also 
suggest an upper bound of a; ~ 0.03 for oxygen clusters to start off. The volume change induced 
by point uranium defects is monotonic but nonlinear, whereas for oxygen defects, increase x always 
reduces the system volume linearly, except dimers that require extra space for accommodation, 
which has been identified as meta-stable ionic molecule. Though oxygen dimers usually occupy 
Willis O sites and mimic a single oxygen in energetics and charge state, they are rare at ambient 
conditions. Its decomposition process and vibrational properties have been studied carefully. To 
obtain a general clustering mechanism in anion-excess fluorites systematically, we also analyze the 
local stabilities of possible basic clustering modes of oxygen defects. The result shows an unified 
way to understand the structure of Willis type and cuboctahedral clusters in U02+i: and /3-U4O9. 
Finally we generalize the point defect model to the independent clusters approximation to include 
clustering effects, the impact on defect populations is discussed. 

PACS numbers: 61.72.J-, 71.15.Nc, 71.27.+a 
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I. INTRODUCTION 



Oxides of the fluorite structure include Zr02, a common ceramic in research and industry, Ce02, and the actinide 
oxides Th02, UO2, and PUO2. The series of actinide dioxides is of great interest in nuclear applications. The present 
generation of nuclear reactors uses UO2 as nuclear fuel. Fast breeder reactors at present employ mixed (U,Pu)02, 
and may in the future use (U, Th)02. In the oxides of the fluorite or CaF2 structure, MO2, each metal atom M is 
surrounded by eight equivalent nearest- neighbor O atoms each of which in turn is surrounded by a tetrahedron of four 
equivalent M atoms. A typical feature of the fluorite structure is the large (|, 5, 5) octahedral holes in which interstitial 
ions can easily be accommodated. Fluorite structure of UO2 transforms to an orthorhombic Prima phase under a 
hydrostatic compression beyond 40 GPa, which in turn followed by an iso-structural transition after 80 GPaji*^ At 
ambient pressure, however, it exists as the single phase, stoichiometric oxide at all temperature up to 2073 K. Above 
that it transforms to the sub-stoichiometric phase U02-a;, whereas at lower temperatures it easily dissolves large 
amounts of interstitial oxygen to form anion-excess compositions U02-i-a;. Higher interstitial concentration leads to 
another ordered phase U4O9, which closely relates to fluorite structure!^ It was argued that stoichiometric U4O9 does 
not exist, and should be U40g_j, actually.'' But for simplicity we still use U4O9 to refer to the non-stoichiometric 
phase hereinafter. There are three polymorphs of U4O9 between room temperature and 1273 K, known as a, [3, 7, with 
the a/ (3 boundary at 353 K and the /3/7 boundary at about 873 K. Only the detailed atomic arrangement in /5-phase 
is clearly determined: the excess anions accommodate in cuboctahedral clusters centered on the 12-fold sites of the 
cubic space group /43d with the uranium sublattice remains undisturhedS^ Although the unit cell is 64 times larger 
than a normal cubic fluorite cell, the average cell is still in fluorite-type except one has to introduce some vacancies 
at normal anions sites and two types of interstitial oxygen, each sited about 1 A from the empty octahedral site of 
the FCC cation sublattice along (110) (O ) and (111) (O ) direction, respectively. This characteristic is also shared 
by the a-phase^ and \]02+xr^ with the difference that U4O9 has a long-range ordering for the interstitial oxygen 
atoms while in UO2+2: it is just short-range ordered. To avoid some oxygens too close together, an intuitive proposal 
that different kinds of oxygen defect are associated to form defect clusters is widely adopted when modeling these 
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phases 7'- 

At first sight the fact that interstitials were detected not at the body centers of the cubic interstitial sites but at 
sites considerably displaced from this symmetric position is puzzling. In rare earth doped alkaline earth fluorides it 
has conclusively shown that at low interstitial concentrations (1 mole % or less) the anions occupy the symmetric body 
center interstitial site, but usually the low-symmetry defect structure is a general feature of anion-excess fluorites:-^' 
About half century has elapsed, people still know few about the stabilization mechanism of Willis O and O sites in 
energetics. In the limit of a; ^ in \J02+x, whether the excess-anions will occupy the octahedral interstitial site or not 
is still unclear. On the other hand, though the occurrence of cuboctahedral clusters in /3-U4O9 has been confirmed by 
experiments, the geometry of defect clusters in low interstitial concentration regime is unknown. One of the simplest 
model is to assume the Willis 2:2:2 cluster (see Ref.Q for its geometry) can exist independently and distribute 
randomly in the material around this concentration. Allen proposed a model for U4O9 in this line by chaining 
2:2:2 clusters along (110) directionii^ Unfortunately his model is definitely wrong because the inconsistencies with 
experimental facts in: (i) leads to an exact stoichiometric U4O9, which might not exist; (ii) no cuboctahedral clusters 
can be formed in his arrangement; (iii) has an equal concentration for O and O sites, against the measurements 
that O position has a much lower occupancy4ii^ 

Therefore to investigate the geometry and stability of possible defect clusters with a first principles method is 
required, but it never be easy. The big unit cell of U4O9 and the shortage of information about atomic arrange- 
ment in UO2+K have restricted most of attempts within point defect approximation, and only formation energy of 
simple intrinsic defects (Frenkel pairs and Schottky defect) were calculated. ^^'^^^yi Applied these energies to point 
defect model (PDM)fi^ii^ however, did not produce satisfactory defect populations — uranium vacancy dominates 
in the hyperstoichiometric regime, against the experimental anticipationj ^^i^^ The failure might be attributed to the 
limitation of the PDM which assumes isolated non-interacting point defects, whereas in UO2+2; this is impossible 
when X > 0.03, as we will show later. Also it can arise from the inaccurate energies produced by the local density 
approximation (LDA) or the generalized gradient approximation (GGA) of the electronic density functional that has 
been proven failed to describe localized states.'^ Nevertheless, some qualitative properties still can access by static 
calculations within this model. For example the diffusion rate of interstitials can be modeled simply by estimating 
the migration energy along all possible paths that bridge the initial and finial interstitial positions, which is readily 
computable by ab initio nudged elastic band (NEB) algorithm. For UO2+X, the conclusion is that a direct diffusion is 
almost prohibited and a normal oxygen on fluorite lattice site must be involved as an intermediate process. That is, 
the interstitial atom pushes a neighboring lattice oxygen into another interstitial site and itself jumps into the vacancy 
thus created {interstitialcy mechanism) The extreme of this process is, evidently, creating a transient oxygen dimer, 
and thus sets up an upper bound to the migration energy for thermodynamical diffusion of oxygens. In order to keep 
the occurrence probability of oxygen dimer being consistent with experimental observation in bulk U409fi^ the energy 
required to form a such kind of dimer should be much larger than the average migration energy. But this has not yet 
been confirmed by ab initio calculations. Near to the surface of UO2 that exposed to air, however, oxygen dimer might 
become prevailing due to oxidations. And their stability in UO2 matrix may shed some light on the mechanism of 
how the material dissolve O2 molecules into individual interstitials. Also, it serves as to verify the Willis' assumption 
that each O interstitial has to be associating with one vacancy that occupies the nearest oxygen sitej since otherwise 
they must form an oxygen dimer. 

These motivate the research work of this paper that mainly focuses on: (i) the stability of isolated point oxygen 
interstitial in UO2+X when x —^ 0] (ii) the stability and decomposition process of oxygen dimer, including the variations 
of energy, cell volume and charges, respectively; (iii) the local stability of defect clusters that composed of oxygen 
vacancies, O and O interstitials. These clusters can be viewed as fractal pieces of a cuboctahedral cluster, which is 
the essential in U4O9 phase. It is believed that the transition from UO24.X to U4O9 involves long-range ordering of 
the defect complexes, leading to a change in the symmetry relating the relative positions of the complexes, without 
producing any atomic re-arrangement within these complexes, i.e., micro-domains of U4O9 should already exist in 
U02+2:.^ What we also want to find out primarily in this paper is what kind of cluster is the most possible candidate 
for this complex, and its polymorphs when x is increased. In next section we will brief the calculation method. Main 
results and discussions are presented in section Hill Sec lIII B] devotes to formation energy analysis and Sec lIII Cl the 
charge transfers, in Sec lIII Dl and IIII El we will discuss the properties of oxygen dimer in UO2 and its decomposition 
process. The defect clustering pattern and its tendency with increased x are given in Sec lIII Fl while in Sec lIII Gl a 
generalization of PDM to include clustering effects is proposed, as well as the associated defect population analysis. 
Finally in Sec lIVI we summarize the paper. 
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II. METHOD OF CALCULATION 

Our investigation on defective behavior of UO2 based on a series of total energy calculations with different configu- 
rations in fluorite structure which varied in simulation cell size and defect arrangement. The plane- wave method using 
density functional theory (DFT) to treat the electronic energy as implemented in VASP code^^^i^i was employed, as well 
as the projector-augmented wave (PAW) pseudopotentials^i^i The 2s^2p^ electrons in oxygen and 6s^6p^5f^6d^7s^ 
in uranium were treated in valence space. The cutoff for kinetic energy of plane waves was set as high as 500 eV to 
eliminate the possible Pulay stress erroneous. Also it has been elevated due to the presence of oxygen which requires 
an energy cutoff at least 400 eV to converge the electronic energy within a few meV. Integrations over reciprocal space 
were performed in the irreducible Brillouin zone with about 8^36 non-equivalent k-points, depending on the system 
size. The energy tolerance for charge self-consistency convergence was set to 1 x 10^^ eV for all calculations. And 
the total convergence of this parameter set was checked well. Without a specific statement, all structures in following 
discussions have been fully relaxed to get all Hellman-Feynman forces (stress) less than 0.01 eV/A. 

The electronic exchange-correlation energy was computed by spin-polarized local density approximation with an 
effective on-site Coulomb interaction to split the partially filled 5/ bands localized on uranium atoms (LSDA-I-U) i^^i^^ 
Parameters of the Hubbard term were taken as J7 = 4.5 eV and J = 0.51 eV, which has been checked carefully for 
fluorite UO2 i^i^i^^i^ Here some comments are desired. It is well known that it is the U but J that contributes to 
electronic structure sensitively. In UO2 case, the value of U quite depends on the atomic arrangement of uranium 
atoms. ^ If uranium sublattice almost being unchanged, as the case here concerned, one can expect the U would not 
vary too much. On the other hand, the influence of interstitial oxygens on localized 5/ electrons should be small 
if they are well separated from uranium atoms. However, as interstitial concentration increased, the impact on U 
may become non-negligible. Therefore we must restrict to certain composition regime, and x < 0.25 should be small 
enough to allow using this set of parameters. This composition value can be estimated roughly by checking the 
induced deformation on the uranium sublattice. The situation of uranium defects is a little embarrassed. Wc cannot 
estimate its effect on U until a more accurate functional becomes generally available, for example, the hybrid density 
functional that has shown impressive versatility in preliminary applicationSi^^ However, for a point defect in a large 
enough cell, to neglect this influence seems reasonable. Another point is about the adoption of LSDA-I-U functional 
instead of GGA+U. The latter has been proven as of a poor description to the defect energetics, which we will discuss 
in details in Sec lIII Gl 

The supercell method has been used to model defect structures. Periodic boundary conditions were imposed on 
the whole system. The geometry of all structures (except those in Sec lIII Fj) are listed in table |T1 where each brick 
indicates a fluorite cubic unit cell (in U4O8), and red points represent oxygen interstitials, which usually occupy 
the cubic centers, except those associated with dimers. Dot-lined box (if drawn) indicates the oxygen cage. No 
atom on the fluorite lattice has been drawn explicitly, except in did and CAidi where the lattice oxygens bonded 
to interstitials were also plotted. Each structure of did, CA2d and CAidi contains one oxygen dimer, respectively. 
Configuration "C81 has the same geometry as C81 but replacing the interstitial oxygen with one uranium, and "C8_i 
or C8_i corresponds to remove one lattice atom from a system with 8 fluorite cubic cells (2x2x2). 

The magnetic effects have been taken into account by initially set up an antiferromagnetic orientation of atomic 
moments. Two cases, the moment ordering along the longest (L) and shortest (S) axis, are considered. The cohesive 
energy Ecoh of each structure is calculated from the total energy by subtracting the isolated spin-polarized atomic 
contributions. Then the oxygen defect formation energy in structure Cm„ is given by 

Ef = Ecoh - rnE^X - '^^02 ■ (1) 

Here m is the number of fluorite cubic cells and n the total oxygen interstitials or vacancies. E02 is the binding energy 
of a neutral dioxygen molecule. Alternatively, one can define an alloy-system like formation energy by choosing Cli as 
one of the reference phases instead of O2 molecule. We call it the relative formation energy, which takes the advantage 
of showing the phase stability of superstructures with different composition explicitly, analogous to that in an alloy 
and compound system i^^i^^ It thus can be calculated as 

Enf-Eco,-{l~^)Eg\-^EZ'^, (2) 

and the value oi n/m stands for the composition of phase Cli in CI, or equivalently, the concentration of oxygen 
interstitial per fluorite cubic cell. All configurations incorporated with uranium defect are marked by a superscript u 
in table m and the formation energy for a defect in "Crrin is defined as 

Ef = Ecoh - mEgl - nEo,u- (3) 

Here Eau is the cohesive energy per atom in the metallic a-U phase, and we use the experimental value of — 5.4eV 
for simplicityi^ 
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TABLE I: Equilibrium properties of uranium dioxide with defects: superscript u denotes uranium defects and negative subscript 
refers to vacancy. AV is the volume difference relative to CI structure and Ef the defect formation energy per point defect. 
Note Ecoh and volume have been averaged to a single fluorite cubic cell. 



Label 



£;coh(eV/cell) 



Volume(AVcell) 



A\/(AVcell) 



E/(eV) 



Structure 



CI 
Cli 

C2i(L/S) 

C4i(L/S) 

C4id(S) 

C42d(S) 

C4idi(S) 

C42(L/S) 



-98.638 
-102.906 
-101.20/-101.199 

-99.71/-99.731 

-99.337 

-100.486 

-100.461 

-101.233/-101.237 



161.34 
157.17 
159.47/159.46 

160.54/160.28 
162.87 
163.05 
162.09 

159.35/159.38 



C8i 


-99.268 


161.05 


C3f(L/S) 


-100.099/-100.361 


160.25/160.16 


C3i'(L/S) 


-101.789/-101.788 


159.34/159.36 


C8-1 


-97.338 


161.54 


''C81 


-98.289 


164.25 


"C8-1 


-96.831 


160.26 



0.0 
-4.17 
-1.87/-1.88 

-0.8/-1.06 

1.53 

1.71 

0.75 

-1.99/-1.96 



-0.29 
-1.09/-1.18 

-2.0/-1.98 
0.20 
2.91 
-1.08 



0.0 

-1.394 
-2.249/-2.248 

-1.413/-1.496 

0.079 

-1.642" 

-1.545" 

-2.316/-2.324 

-2.169 
-1.509/-2.294 

-1.853/-1.850 
7.525 
8.194 
9.056 
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" per two oxygen interstitials 



Vibrational frequencies of interstitial oxygens were calculated by finite difference method with frozen phonon ap- 
proximation. At finite temperatures, these vibrational frequencies contribute to the first order of defect free energy 
directly, which is given by F{T) — Ef — kbTIh Zy, with the partition function 



3=0 



exp 



kbT 



(4) 



where is the Boltzmann constant and Ej the eigenvalue energy for the j-th vibrational mode with frequency w*, 
and the harmonic approximation Ej = fujj'^{j -f i) has been used. Mind here we have not subtracted the vibrational 
free energy of the reference state O2 molecule, and comparison of the calculated free energies therefore can be made 
only among configurations with the same number of interstitial oxygens. 

Regarding charge transfer calculations, it is well known that the concept of static atomic charge in ah initio calcula- 
tions usually leads to ambiguity due to the arbitrariness in determining the belongingness of electrons. Nevertheless, 
there are several methods exist to compute the effective atomic charge, which do provide some useful qualitative 
understanding. Among those the Bader's conception that to partition an electronic density by surfaces formed by the 
density minimums (zero flux surfaces) is one of the most intuitive. It is simple to calculate Bader charges, requiring 
only atomic positions and electronic density as input. The partition surfaces are determined by finding the charge 
density minimums.— Then the atomic charge is obtained by subtracting the valent electrons from the integral of 
charge density over the space surrounded by the partition surfaces that envelops the atom. Another widely used 
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FIG. 1: (Color online) Calculated variation of the volume change in UO2+1 with the deviation from stoichiometry x: Solid 
squares stand for point oxygen defects, and the solid line is the linear fitting to them; open triangles are for uranium defects; 
solid circles are those incorporated with one oxygen dimer; the dotted line is for eye guide. Experimental data: dashed line is 
for U02-i-a; reported by Alekseyev et al. and others for /3-U4O9 (at room temperature M, 503 K ffl, and 773 K ©, respectively.) 

concept is dynamical effective charge, defined by the change of polarization induced by atomic displacements,"^^ which 
is beyond the scope of this paper and will not elaborate here. 

III. RESULTS AND DISCUSSIONS 
A. Dioxygen molecule 

We first discuss the dioxygen molecule. The O2 molecule was modeled by putting it in a periodic cubic cell with a 
lattice constant of 15 A, large enough to eliminate the factitious interaction among its images. Only one k point (F) 
was used. Since the notorious failure of LDA in description small isolated molecules, we employed here (and only here) 
the revised Perdew-Burke-Ernzerhof (rPBE)''^ GGA electronic exchange-correlation functional. The bond length was 
optimized to be 1.22 A, in a good agreement with experimental 1.21 Ai^i^ The calculated binding energy is — 5.75eV, 
a little deeper than observed — 5.1eV<^ This discrepancy should attribute to the difficulty of current functional to 
take into account the van der Waals interactions accurately. The vibrational frequency of stretch mode, however, 
was well reproduced as 1588. 6cm~-'^, against the experimental 1580.2cm~^i^ As a check to the validity of Bader's 
conception, we calculated the Bader atomic charge for each oxygen atom in O2 and got them as ±0.09 e, reflects the 
essential of covalent bond correctly. The deviation can be reduced further when in an ionic bond environment where 
the charge density minimum surfaces sharply show up. 

B. Structure and formation energies 

1. Oxygen interstitials 

The calculated equilibrium properties of 14 configurations, including cohesive energies, equilibrium volumes, volume 
changes relative to the ideal UO2 cell and defect formation energies are listed in tableU These data have been averaged 
to one fluorite cubic cell. It can be seen that the cohesive energy always decreases as oxygen interstitial concentration 
increased, demonstrating the tendency of uranium dioxide to dissolve oxygens. The solubility, however, cannot be 
determined by simply taking the limit of this cohesive energy vs concentration curve. Also, the relative stability 
among different configurations has been obscured here. To get that information explicitly, one needs back to the 
relative formation energy. 

One interesting thing is that we flnd except that of oxygen dimer, introducing point oxygen interstitials always 
shrinks the system, i.e., leads to a negative AV, as shown in figure [T] This feature differs from GGA resultSfi^ 
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FIG. 2: (Color online) Formation energy of oxygen interstitials in UO2 arranged in various configurations. C81 corresponds to 
an isolated defect approximation and other configurations must be interpreted as ordered defect phases. 



but agrees with GGA+U;^ and may attribute to the behavior of localized 5/ electrons. Generally, a negative AV 
means the interaction between the matrix and the interstitials is dominated by attractive chemical potentials rather 
than by mechanical effects (atomic size effect). The latter always results in a swollen volume and is important for 
big interstitial atoms or inert gases. Oxygen dimer belongs to this class and requires extra space to accommodate, 
which can be seen more clearly when compare C4if; with C4i and C\\d\ configurations. The influence of magnetic 
orientation on equilibrium volume is almost negligible except in the cases of C4i and CSf, of which only C3f has a 
notable formation energy difference between L and S orientation. 

The calculated slope of volume variation induced by oxygen interstitials (the solid line in figure [1]) is in a good 
agreement with experimental change of the lattice constant a = 5.4696 — 0.1495a; as reported by Alekseyev et al. for 
homogenous UO2+2; powders as quoted in Ref. 'ss'l (the dashed line). Also, it is in accord with the volume change of 
/3-U4O9 measured at room temperature^^ with respect to that of stoichiometric U02-^ Increase temperature to 503 K 
and 773 K expands the material greatly,-^ which can be understood in terms of thermal vibration effects and extensive 
defects generation. 

Figure [2] shows the defect formation energies of oxygen interstitial in all considered configurations of UO2+2: within 
< a; < 0.25. Note the value of C42d and CAidi are for two interstitials. A remarkable feature in this graph is that in 
energetics an oxygen dimer mimics a single oxygen atom. Comparing that of the perfect crystal CI with C4id, and 
that of C4i with C42d and CAidi , we see that despite the latter contains one more interstitial, the formation energy 
is almost the same. That means to absorb an oxygen from O2 gas into UO2 and forms a dimer will neither release 
nor gain heat. Point interstitial and dimer would have almost the same behavior except that a dimer needs a bigger 
space for accommodation. This mimic is also supported by Bader effective charge calculations: they have almost the 
same charge too (see below). However, this does not suggest the stability of oxygen dimers in UO2 since point oxygen 
interstitial always has a lower per atom formation energy. 

Our calculations also present a remarkable system size dependence in formation energy, contrasts to that of GGA 
results where a value of —2.6 and — 2.5eV were obtained for Cli and C2i configurations (almost size- independent), 
respectively)^ revealing the limitation of applying the pure GGA to defects in spite of its impressive performance in en- 
ergetics of perfect bulk U02.^ No magnetic ordering and volume relaxation were considered in that GGA calculation 
A discrepancy about 1.2 eV with our result for Cli , however, cannot be covered by these effects since volume relaxation 
would definitely increase the discrepancy and magnetic contribution can not be of that magnitude, and it therefore 
should attribute to the behavior of localized 5/ states. 

The deepest formation energy shown in figure [5] is — 2.32eV (configuration C42), rather than the isolated approx- 
imation of a point interstitial's (C81) — 2.17eV. Actually, except those configurations with 8 fluorite cubic cells, the 
defects in all other structures cannot be interpreted as isolated ones because the non-negligible interactions among 
their images arisen from periodic conditions. This invalidates the defect stability analysis based on their formation 
energy directly. Mapping these configurations onto an alloy system can circumvent this difficulty. Namely, to view 
these configurations (discard those with dimer) as an alloy system with oxygen interstitials distributing over the 
fluorite cubic centers (U408-I-;/). Then the extreme phases of this system are CI and Cli. Following this way, figure 
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FIG. 3: (Color online) Relative formation energy of different phases in U40g+y. C8i corresponds to an isolated defect 
approximation and other configurations must be interpreted as ordered defect phases. 



[2] transforms into figure [3] with the help of Eq. ((21) , where the sold fine indicates the ground state hull. We then 
find that C3f (S) and CSj are close to be ground states, while C8i, the isolated point interstitial approximation, 
will decompose into a mixture of CI and C42 phases. That means, defect clustering is inevitable when x > 0.03.— 
Since C42 may not be the physical ground state (mind neutron diffraction experiments suggested that no octahedral 
site should be occupied around this compositionii^,) , this limit can be lowered further. On the other hand, it seems 
reasonable to assume that C8i has already approached the limit of an isolated point interstitial, namely, no notable 
formation energy would be gained or lost if enlarge the system to C27i or C64i. If it is true, then the isolated point 
interstitial will always exist when x — > 0. Its site, according to structure symmetry, should be the octahedral position. 
It is worthy to note that the PDM fails at about x ~ 10^** instead of 10^^ with GGA formation energiesi^^ This two 
orders discrepancy is owing to the inaccuracy of the formation energies they used, which can be improved greatly by 
LSDA+U method, see Sec lIIIGI for details. 



2. Other defects 



Point oxygen vacancy and uranium defects are all modeled in a system with 8 fluorite cubic cells. Namely by C8_i 
and "C8±i. The volume change in C8_i is in accord with that of point oxygen interstitials: linearly decreasing with 
an increased x, fitting to /\V = 0.01 — 14. 7x, as shown by the solid line in figure [T] Uranium vacancy also obeys 
this law, whereas the interstitial has a much rapid change. Totally, they still decrease monotonically with x, but no 
longer in linear. All three defects have a formation energy larger than 7eV, in contrast to previous ab initio results. 
We prefer to postpone this discussion to Sec lIII Gl with Frenkel and Schottky defects together. 



C. Charge transfers 



It has long been believed that dissolve oxygen in UO2 will oxidize U"'+ to U'''+, even state. The exact charge 
transfer induced by defects, however, is unclear. Qualitative analysis is accessible to this problem with empirical 
shell modelf^ nevertheless the calculated energy depends on atomic positions sensitively,^^ obscuring its applicability 
to defects with noticeable structure deformations. A direct calculation of the charge state from first principles is 
therefore desired. 



1. Oxygen interstitials 



The calculated Bader effective charges using electronic density generated with VASP code are listed in table [TTl 
where the interstitials and the lattice oxygens that forming a dimer are excluded from the average operations, and 
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TABLE 11: 


Bader effective charj 


res of UO2 with defects: average 


charge q, standard deviation a, difference from that 


in perfect 


UO2 Sq and the maximal transferred charge A; 


(±0.02) 


. All 


in the unit of posit: 


ron charge 


e. 






Label 


defects 




uranium 








oxygen 






q 


q 


O" 


Sq 


A 


9 


(T 


oq 


^ TJi a X 


CI 




2.56 


0.0 


0.0 


u.u 


1 Oft 
-l.Zo 


u.u 


u.u 


0.0 


Cli 


-1.04 


2.62 


0.11 


0.07 




lift 
-1.10 




10 

U. lU 


0.11 


C2i(L/S) 


-1.15 


2.63 


0.11 


0.08 


u.z / 


1 0/1 


n no 
u.uz 


HQ 
U.Uo 


0.07 


C4i(L/S) 


-1.18/-1.14 


2.60 


0.08 


0.04 


u.zo 


1 Ofi 
- 1 .ZD 


U.Uo 


09 
U.UZ 


0.08 


C4id(S) 


-O.e^-O.??)" 


2.56 


0.02 


0.01 


-U.Uo 


1 Oft 
-l.Zo 


U.Ul 


n 
u.u 


-0.03 


C42d(S) 


-0.66,-0.59 
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listed separately in the defects column. We find all oxygen interstitials that occupy the cubic body center having 
a charge state close to the lattice oxygens, especially in the C81 phase where the difference is only 0.03 e. In C81 
the disturbance to lattice oxygens is also small, the largest charge transfer is just 0.05 e. A similar situation holds 
for uranium atoms, except two of them lost about 0.24 e, respectively, which contribute to the standard deviation 
directly. Considering oxygen and uranium in perfect UO2 have only a charge of —1.28 and 2.56 e, all are smaller 
than the nominal chemical valences but close to that of a partially ionic model that widely used in semi-empirical 
potentials ,'^"' we can re-interpret the Bader charges by multiplying a scaling factor to make them comparable with 
the chemical valences. In this sense the change of the charge state in these two uraniums should be about 0.5 e, i.e., 
they are oxidized to U'*'^+ instead of U^+. The transferred charges, however, cannot cover the amount absorbed by 
the interstitial oxygen, and all other normal uraniums and oxygens have also lost a small portion of their charge. 
This observation contrasts to the conventional expectation and reveals the difficulty to oxidize uranium to a higher 
valence state. The charge transfers in other configurations also support this point: in all cases each oxygen interstitial 
can oxidize two and only two uraniums to U'*'^^ while leaving others almost unchanged, no higher valence state of 
uranium has been observed. As to which uranium is apt to be oxidized, obviously the answer is the nearest neighbors 
(NN) of the defect, but oxidization of some next NNs also was observed. It is worthwhile to point out that we did 
not find a sensitive dependence of the charge state on Hubbard U parameter. 

The more deformed the geometry is, or equivalently, the more interstitials the system contains, the charge state of 
lattice atoms are disturbed more drastically. It is clear by comparing the charge transfers in C3f with CSj , or C4i 
with C4:2. The largest Amax for oxygen takes place in Cli with the largest composition, and in C42d with a dimer. 
The smallest A.max for uranium and oxygen are in CAid, also containing a dimer, both are — 0.03e. The difference 
between CAid and C42d is the former contains only one interstitial which bonds to a lattice oxygen and the latter 
contains two interstitials that bonding to each other. Table HIl illustrates that in the former case no charge has been 
transferred from other lattice atoms, and only charge redistribution within the dimer is involved that making it has 
a total charge close to a lattice oxygen; in the latter case, however, absorbing charges from other atoms is necessary 
and gives them a similar charge state as the interstitial in C4i, especially only two uraniums are oxidized to U'*'^''' 
state here in spite of there are two interstitials presented. The total charge of the dimer, — 1.25e, close to a lattice 
oxygen in UO2, indicates it should be Ol~ actually. 

It is worthwhile to note that oxygen changes its charge state almost continuously but it is discrete for uranium 
when lost its charge. That is, except those atoms who lost ~ 0.25 e, the changes of charge in other uraniums are less 
than 0.03 e. Moreover the discrete lose of charge is always accompanied by lowering the local moment of uranium from 

2 /iB to 1 /iB • Since the local moment of uranium in UO2 originates from localized 5/ states, it is obvious that 5/ 
electrons contribute to this process greatly. This can be understood in the partially ionic charge model: although the 
chemical valence of uranium in UO2 is 4-l~, table HIl shows in fact the physical valence has only 2.56-f . Namely, only 
the 7s^ and a fraction of 6d^ electrons are completely transferred to oxygen. Uranium cation still holds about 0.24 e 
of the 6d^ electron and other remainder forms two weak U-0 covalent bonds, each has a portion of ~ 0.2 e. When 
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FIG. 4: (Color online) The diflterence charge density of an oxygen dimer in C42d configuration projected onto the [100] plane 
crossing the dimer center. An analogous density holds for the dimer in C4i£j. 

oxidized by oxygen interstitials, the cation loses its 6d^ electron completely (transferred to the interstitial atoms). As 
a consequence one of the localized 5/ states becomes the outermost orbit, which spreads extensively and eventually 
the cation lost half of its local moment. This mechanism also explains the difficulty to oxidize uranium to a higher 
charge state since transferring a 5/ electron requires much larger energy than 6rf one. 

2. Other defects 

In point oxygen vacancy case (configuration C8_i), the uranium cations gain charges and decrease the average 
valence to 2.53+, but the disturbance to remain oxygen is small. The largest charge transfer for uranium is — 0.34e, 
associating with three other uranium atoms each of them gets an extra charge about —0.25 e, respectively. Compared 
with the interstitial case, here no notable change in local moment was observed. The value of —0.34 e implies the 
cation has retracted the portion of electrons shared by the removed oxygen (~ 0.1 e), and the — 0.25e indicates that 
each quarter-filled states of 6d electrons seems stabler than continuous occupancy. 

Point uranium vacancy is analogous to two oxygen interstitials in that there are four uranium cations lost the charge, 
three NNs and one next next NN, ranging from 0.23 to 0.26 e. All of them also lost half of their local moments. The 
change in other uraniums is negligible. However, it disturbs the oxygens severely, with a A^ax as high as 0.13 e, even 
though the averaged charge is still close to the perfect one. The oxygen charge state in "C8i is almost the same as in 
C8_i, except that here there are six (NNs) instead four uraniums gain charge, ranging from —0.19 ~ — 0.25e. Again, 
no apparent impact on other atoms. The extra charge provided by the interstitial uranium is almost absorbed by its 
six NNs completely. The magnetic ordering has been damaged severely, and the change in exchange interaction has 
made some 5/ electrons flip their spins, but no uranium was observed to have a moment of ~ 1 /is- 

D. Oxygen dimer in UO2 

As previous sections mentioned, although oxygen dimer has a similar behavior in energetics and charge state as a 
single oxygen interstitial, it actually is an ionic molecule, and formed when oxygens are forced to close to each other 
enough. But this is difficult due to the energy barrier between the individual atoms. In UO2, irradiation provides 
enough excess energy to overcome this barrier. For example in an a decay the recoil of the daughter nucleus produces 
a ballistic shock with an energy release of about 70 keV^^i^ which frequently takes place in nuclear fuels. Nonetheless 
this cannot survive the dimers to equilibrium conditions, even if they do appear transiently. Another situation where 
oxygen dimers can be observed is near the surfaces exposed to oxygen gas. Oxygen molecules adsorbed onto the 
UO2 surface will obtain additional charges then diffuse inwards. Decomposing the molecule at the vacuum side of 
the surface is almost impossible due to the large binding energy, while in UO2 side it prefers to oxygen sublattice 
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FIG. 5: (Color online) Behavior of oxygens during a dimer decomposing process along (111): potential shape (solid lines) and 
equilibrium intra-distance (dashed line), 5d is the initial separate distance and Ad the final bond-length; inset: changes in cell 
volume (black line) and atomic charges (red lines). 



sites instead of the interstitial positions, where it decomposes into individual interstitials, with a barrier only about 
0.21 eV (see below) 4^ 

Figure [¥] shows the difference charge density (reference to the corresponding atomic charge) of an oxygen dimer 
in UO2 (C42d) projected to [100] plane, as the arrow indicates. The covalent bond between two interstitial oxygens 
presents evidently. A similar picture has been observed in C^id configuration or a natural O2 molecule. Analysis 
shows that it in fact is an 0\~ , with the two additional electrons occupying the 2p7r* antibonding orbitals and the final 
bond order is one. The calculated bond length is 1.39 A, slightly short than experimental 1.49A.'^'^ This discrepancy 
owes to the compression from oxygen cage and can be removed. For example when the dimer is formed by bonding 
to one lattice oxygen (C4irf), where although the charge state is still similar (— 1.38e), the bond length extends to 
1.47 A, in a good agreement with experimental data. 

Accommodation oxygen dimer in UO2 leads to a swollen of the system volume (see table The induced stress 
forces them orientated in (111) direction and occupy the Willis O sites actually. In energetics, oxygen dimer in UO2 
is meta-stable, see figure [51 Its decomposition process can be modeled by successively moving the interstitial oxygen 
(as a test atom) in C4id along the (111) direction till to the cubic center, which is the most possible separate path. 
The resulted potential shape is shown in figure [H where (5d is the initial depart distance between the two oxygens 
and Ad the final (dimer) length. The structure frozen line was obtained by fixing the cell and all other atomic 
positions whereas the optimized one resulted from a fully relaxation of the cell volume and shape and the nearby 
atomic coordinates that surrounding the defect. 

Note a distance of 5d = 2.2 A stands for the state that the initial position of the test oxygen already close to the 
cubic center. From figure[5]we get the critical distance to break a dimer is about 1.73 A, with a barrier of 0.21 eV. Inset 
gives the variations of system volume and Bader atomic charges of the two oxygens, demonstrating a drastic behavior 
around the breaking point. Two points need to be noticed here: the large charge transfers and the contraction of 
system volume. The later confirms that atomic size effect is not an important factor for oxygen interstitials in UO2 
where chemical interaction overwhelming. A deduction of this is that a single oxygen interstitial can occupy a site 
other than the cubic center, regardless it has the largest space. Indeed, no experiment has detected the occupation 
of this site in U02+:!; when x > 0.1. Chemical interaction might prefer other sites if volume is expanded. As the 
"structure frozen" line shows, interstitial oxygen is apt to forming dimers when the volume is fixed at 651.49 A . 
Therefore, oxygen dimers may also exist at regions with negative stresses. 

As the limit case of an interstitialcy diffusion model, to form an oxygen dimer in UO2 requires an energy of ~ 1.75 eV 
(figure [5]), compatible with the NEB migration energy of 1.1 eVJ^ This magnitude of migration energy corresponds 
to 6d ~ 1.8 A, with an equilibrium intra-atomic distance about 2.0 A and the atomic charges — l.Oe. Therefore a 
charge oscillation induced by oxygen diffusion is about 0.2 e, almost the same level as to oxidize one uranium. 
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TABLE III: First principles results for structural, energetic, and vibrational properties of oxygen interstitial and O2 dimer 
in different configurations. For comparison, calculated values for O2 in vacuum are also listed. AE is the energy difference 
between interstitial O/O2 and vacuum O2 (formation energy per pair interstitials), do is the equilibrium bond length, q is the 
Bader effective charge, uj is the harmonic frequency. Note the q in the last row is just for to label the experimental condition. 



Label 



AS (eV) 



do (A) 
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-4.496 
-2.993 
0.159 
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1588.6 

1580.2V- 



" Ref. 
Ref. 



E. Vibrational frequencies 



Raman and infrared spectroscopies provide information about atomic vibrations. These techniques can be em- 
ployed to detect defect clusters by searching the characteristic vibrational frequencies. At finite temperatures, these 
frequencies directly contribute to the formation energy and structural thermodynamic stability. 

Vibrational frequencies of single oxygen interstitial (has three modes) and dimer (has six modes) in C2i, C4i, 
C4,id and C4,2d configurations were calculated. In all calculations, we aligned the magnetic ordering direction along 
the shortest axis (S), which always has the lowest energy. Only harmonic frequencies were computed here and has 
omitted all anharmonic effects. For fluorite UO2, we have checked that the contribution from the latter is very small 
for oxygen and uranium interstitials (within ±3cm~^). Table Hill lists the calculated frequencies (w), as well as the 
equilibrium bond length for dimers (dg) and formation energies {/S.E). Due to the compression from oxygen cage, the 
vibrational frequencies in C4,2d have greater value than their counterparts in C4id. The stretch model of O2 molecule 
(with the largest w) has been greatly softened when incorporated in UO2. This is analogous to the incorporation of 
H2 in an interstitial position of semiconductors)^ where a decrease of the binding energy, an increase in the bond 
length, and a lowering of the vibrational frequency was observed. The underlying physics, however, might be different. 
In this case, by comparing the calculated Bader effective charges with the partially ionic model of U02,'*^ we can 
identify the nominal charge of the oxygen dimers should be about —2.0 e. The variation of bond length confirms this 
interpretation. Consequently, the frequency of stretch model is lowered from 1588.6 cm^^ to 995.4 cm~^ in C4:2d and 
795.4 cm-i in CAu- 

According to the calculated static energies, C42d will decay to C2i, and CAid to C4i eventually (see /\E in tablc lllip . 
Computed frequencies indicate thermal vibrations would accelerate this process further. Figure [5] gives the difference 
of free energy between CA2d (Fd), C2i {Fg), CAid (fd) and C4i (fg) calculated with their formation energies (table 
H]) and vibrational frequencies (table IIIip , respectively. The rapid drop of the free energy differences with increased 
temperature implies that meta-stable oxygen dimers in UO2 have a very short lifetime at finite temperatures, and 
with little possibility to occupy the cubic center sites: they must have been decomposed before enter the oxygen cage. 



F. Defect clustering in UO2-1-; 



This section devotes to the possible defect clustering pattern in U02+a;- Instead of compute the formation energies 
directly, we focus on the local stability of O and O sites in different circumstances here. This method cannot 
determine what cluster is the most stable one, but it does rule out some combinations of O , O and oxygen vacancies. 

For this purpose we calculated the potential landscape felt by a test oxygen atom. Just one fiuorite cubic unit cell 
was used. Here since the cell and all atomic coordinates have been frozen up except that of the test oxygen, the error 
introduced by periodic conditions is proportion to the second order of the charge density variation 6p that induced by 
images of the test atom. This precision is enough for a qualitative analysis (Mind ionic interactions among the test 
atom's images contribute only a constant to the energy and therefore irrelevant to the problem). 
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FIG. 6: Variation of free energy difference contributed from interstitial vibrations. Inset: the free energies of a dimer in C42d 
and its relative stable state C2i as a function of temperatures. 
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FIG. 7: (Color online) Potentials for an oxygen interstitial in UO2 along (111) (left) and (110) (right) direction crossing the 
cubic center. The numerics refer to lattice constant and arrows point to the position of O or O site. Fractional coordinate 
0.5 denotes the cubic center. 



1. Local stability of basic clustering modes 

At first we check the local stability of a single O and O site. Figure [7] shows the potential shapes crossing these 
two sites. It is seen that O site becomes meta-stable when the lattice constant increases to about a = 5.44 A. And 
isotropic expansion stabilizes this site further which makes it the global minimum if a > 5.6 A. A single O interstitial 
actually forms a dimer with the nearest lattice oxygen and this behavior is in consistence with the structure frozen 
curve in figure O However, this effect does not benefit the stabilization of O site. Under ambient conditions the 
experimental lattice constant for UO2+X is within 5.45 ~ 5.47 A, therefore a single O or O oxygen interstitial (as 
well as clusters formed by them only) is almost unstable. 

The simplest cluster involved one oxygen vacancy, say a V-0 or V-0 pair, is obviously unstable since nothing 
can prevent them from annihilation. The next triple cluster is V-0 (O ) pair stabilized by an O or O interstitial. 
Considering the short distance between the nearest O and O sites, the situation should be quite similar for them. 
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FIG. 8: (Color online) Potentials for an oxygen interstitial in UO2 (incorporated with a V-0 pair) along (111) (left) and (110) 
(right) direction crossing the cubic center (fractional coordinate 0.5). 
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FIG. 9: (Color online) Potentials for an oxygen interstitial in UO2 (incorporated with two V-O pairs) along (111) (left) and 
(110) (right) direction crossing the cubic center (fractional coordinate 0.5). 



Therefore hereinafter we only consider the cases that incorporated with V-0 pairs. The potential shape for an 
O (O )-V-0 cluster was calculated in an analogous manner except a lattice oxygen (0.75, 0.75, 0.75) was moved to 
(0.883,0.5,0.883), a nearest O site, to create the V-0 pair, as shown in figure [8l Although the curve along (110) 
already changes asymmetrically about the cubic center (with a fractional coordinate 0.5), O is locally unstable since 
it will decay to O (with a swallow trap) then finally to a position beyond the (0.75, 0.75, 0.75) site. This rules out the 
O -V-0 (V-20 ) and O -V-0 (V-20 ) triple clusters that distribute symmetrically about a lateral of the oxygen 
cage. 

To stabilize the O site locally, we have tried all possible combinations and find two nearest oxygen vacancies seems 
necessary. Figure [3] gives the potentials that incorporated with two V-0 pairs. Mind these O sites should be in 
otherwise empty oxygen cubes that do not share the lateral linking the two vacancies with the original one. The pairs 
arc thus created by moving (0.75,0.75,0.75) oxygen to (0.617,1.0,0.883) and (0.75,0.75,0.25) to (1.0,0.617,0.117), 
respectively. We see it do prefer the O but not O site. In fact this cluster would become the Willis 1:2:2 (O :V:0 ) 
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FIG. 10: (Color online) A schematic phase diagram for oxygen defect clustering in \J02+x- The boundaries, however, are not 
clearly known. 

cluster— if move the two O interstitials to their nearest O sites and form two V-0 pairs rather than the V-0 ones. 
Figure [9] shows it might be local stable, in consistent with empirical calculations.— The stabilization of O by V-0 
pairs sited in the otherwise empty oxygen cages is unclear by figure[8l but calculations show a V-20 triple do stabilize 
O" locally (0"-V-20'), as well as a V-20" triple can (0"-V-20"). 

Thus we finally arrive at the conclusions that: (i) O or O interstitials cannot exist by themselves; (ii) each O 
site must incorporate with two nearest oxygen vacancies, while O can be stabilized by a V-20 (O ) triple. That 
means the possible clustering pattern for oxygen defects should only be: (a) V-30 or V-40 isolated clusters, in 
the same manner of split-interstitial where several atoms sharing a single lattice site; (&) cluster chains of V-0 -V or 
V-20 -V by sharing the vacancy sites. These chains should be closed to have all O interstitials locally stable whilst 
minimize the vacancy/interstitial ratio; (c) cluster of V-(2)0 -V chains terminated by two V-(2)0 clusters at both 
of the extreme sides by sharing the vacancy sites. We call these small fractal clusters the Willis type cluster, including 
1:2:2, 2:2:2,**-^" 4:3:2^*'-'*^ clusters and so on. However, their actual stability is still unknown which requires accurate 
knowledge about their formation free energies. 

2. Phase diagram for clusters 

In U02+a; , the positive formation energy of oxygen Frenkel pair and the small energy gain from interactions among 
interstitials (see figure [2]) implies the only way to reduce the energy increment from creating vacancies is via vacancy- 
interstitial (V-I) interactions. Obviously the nearest V-I pairs contribute the most. Therefore the relative stability 
of clusters can be judged roughly by counting the number of nearest interstitials around each vacancy. For example 
in a 1:2:2 cluster, each vacancy has only two V-I pairs, while in 2:2:2 it has three and in 4:3:2 it is 3.3. That means 
1:2:2 should be meta-stable, even though it can explain the concentrations measured by Willis in 1964?^ However, this 
data also can be explained by a larger cluster with 4 O interstitials, namely, a 2:2:4 cluster where each vacancy has 
four V-I pairs. Willis type clusters are necessary in order to explain the large concentration of O interstitials, which 
is impossible by cuboctahedral cluster only (belonging to pattern (&)). For example the data for crystal A done by 
Murray et al.*^ obviously belongs to 2:2:2 clusters while crystal B should be a mixture of 4:3:2 and cuboctahedral 
clusters or a 6:4:2 cluster. 

However, a big Willis type cluster is unfavorable since the disturbance to fluorite lattice is proportion to its size 
linearly. A similar situation holds for a loosely closed chain of the pattern (&). In this sense the most regular and 
close-packed defect cluster, the cuboctahedral cluster, takes the advantage of sharing the space with all vacancies and 
interstitials to minimize the damage to the matrix. Also, one fluorite cubic cell can accommodate one (or less) Willis 
type cluster or one cuboctahedral cluster. But the former provides only 2 excess anions while the latter provides 5. 
When composition x increased, cuboctahedral cluster has a big advantage over Willis type cluster, not to mention 
each of its vacancy has a number of V-I pairs greater than three. As for the clustering pattern (a), though there are 
3 (V-30 ) or 4 (V-40 ) V-I pairs for each vacancy, we can discard them at first since no experiments showed so high 
concentration for O interstitials. 

It becomes evident when check the variation of O :0 ratio as x increased: around a; = 0.11 ~ 0.13, three different 
data sets were observed (0.08:0.16,^ 0.14:0.12, and 0.33:0. 10^^), implying the occurrence of Willis type clusters. As x 
approaches to 0.25, however, this ratio increases drastically^!^ and shows the dominance of cuboctahedral clusters. 
Therefore, by taking the stability of point interstitial at low x into account, one concludes that there is a quasi- phase 
diagram for oxygen clusters in U02+£c , as shown in figure [TUl^ To determine the exact geometry of ground clusters 
and their boundaries would be the center of future works in this field. 

It is worthwhile to point out that such kind of defect clustering is not unique to uranium dioxide. According to the 
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formation energy of point defects, one can classify binary compounds into three classes: (A) all formation energies 
are positive, (B) only one of the formation energies is negative, and (C) both cation interstitial (vacancy) and anion 
vacancy (interstitial) have negative formation energies. There is no off-stoichiometry driven force in case (A), and 
disfavors extensive defect clustering. However the negative formation energies in the other two cases will drive the 
system to a non-stoichiometric composition where defect clustering becomes favorable. This is because the interaction 
among defects can lower the system energy greatly, and lead to a pure defect clustering (via a full vacancy-interstitial 
annihilation) or mixed defect clusters that contain both vacancy and interstitial. Also, the mixed cluster is possible 
only when the point defect with positive formation energy (vacancy or interstitial) has the function to stabilize the 
other defects in an energy favorable configuration (in a similar concept of the split-interstitial defect mechanism). 
Obviously UO2 fulfills these conditions (see figure [2] and previous subsection) and belongs to case (5), where oxygen 
interstitial has a negative formation energy and clustering involves no uranium sites. On the other hand, case (C) 
contains clusters composed of both cation and anion defects, and might exhibit more complex behaviors. 



G. Concentration of defects 

1. Generalization of the PDM 

The point defect model (PDM) was introduced by Matzke^^ and Lidiard^^ to analyze the populations of defects in 
UO2+2;, where x indicates the deviation from stoichiometry. This model is based on the hypothesis that the defects 
responsible for the deviation from stoichiometry in UO2+X sue isolated point defects. However, it has been known 
for long that oxygen interstitials form clusters, and usually PDM performs poorly even at small |a;|.— Therefore 
it is worthwhile to generalize this model beyond the point approximation. Since traditionally defect concentrations 
are defined in a lattice model as the number of defects present divided by the number of available sites for the defect 
under consideration, the most general and elegant generalization of PDM would be cluster variation method (CVM),— 
which also bases on lattice gas model and computes cluster configurational entropy explicitly. The related effective 
cluster interactions can be determined by cluster expansion method (CEM)i^ For UO2+X, at first sight it seems being 
a quaternary system (Vb, Vu, lo, and lu) and cannot be tackled by modern CVM and GEM techniques. But since 
defects on the uranium subsystem usually are isolated point defects that couple with oxygen subsystem via Schottky 
defects, in fact only oxygen defects need to be treated explicitly. However, in order to include O and O sites in the 
calculation, one has to use an extended lattice, which introduces another two difficulties. 

The first one relates to the local stability of O and O sites, since these sites are not well defined and usually a full 
relaxation is required to get the optimized geometry. However, in most configurations they are not at the potential 
minimums and makes it impossible to include the relaxation effects in the ab initio GEM procedure. Fortunately, an 
algorithm proposed by Geng et al. can tackle this problem simply.^ - The second difficulty is that most configurations 
on the extended oxygen sublattice is unphysical, i.e., some distances among oxygen sites are too short to allowed. To 
exclude these unphysical configurations, one has to use loose clusters to expand the energy, which deteriorates the 
convergence of cluster expansion drastically. 

If all non-negligible clusters are independent and uncorrelated, a simple approximation exists to calculate cluster 
populations. Two clusters are called independent if none of them is the other one's sub-cluster (or loosely, cannot 
dissociate or combine into other clusters). This ensures all cluster concentrations are completely independent. Assume 
there are M such kind of clusters under consideration, then the system free energy can be written as 

M 

F = Y,p,{E, + KBT\np,) (5) 

i=l 

in the closed regime (in which the system cannot exchange atoms with the exterior). Here Ei stands for the z-th 
cluster's formation energy. Minimization this free energy with respect to each cluster density pi (under the condition 
that X is fixed) gives 



associated with the composition equation 



Pi = gi exp ( — ^ ) , (6) 



/(Pi, • • ■ ,Pm)- (7) 



In Eq.([5]) the factor gi is introduced to account for the degeneracy if the cluster has internal freedom, while non- 
degenerated states can be treated as independent. This gives the internal entropy contributions and is the most 
significant difference between the independent clusters approximation (ICA) and PDM. 
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TABLE IV: Formation energies (eV) of point defects in UO2: uranium and oxygen vacancies (U-Vac and 0-Vac), uranium and 
oxygen intcrstitials (U-Int and 0-Int), Frenkel pairs (O-FP and U-FP), and Schottky defect (S). 





U-Vac 


0-Vac 


U-Int 


0-Int 


O-FP 


U-FP 


S 


LSDA+U 


9.1 


7.5 


8.2 


—2.2 


5.4 


17.2 


10.6 


GGA+U'' 


8.4 


4.5 


4.7 


-0.4 


4.0 


13.1 


5.8 


GGA'^ 


4.8 


6.1 


7.0 


-2.5 


3.6 


11.8 


5.6 


GGA"^ 


5.1 


6.1 


7.5 


-2.6 


3.5 


12.6 


6.0 


LDA" 


3.3 


6.7 


7.3 


-2.9 


3.9 


10.7 


5.8 


LDA-LMTOii 


19.1 


10.0 


11.5 


-3.3 


6.7 


30.6 


17.1 


Semi-empirical^ 


80.2 


16.9 


-60.8 


-12.1 


4.8 


19.4 


11.3 


PDM estimateai^ 










3.0-5.8 


9.5 


6.0~7.0 



" this work, with 8 fiuorite cubic cells 
^ with 8 fiuorite cubic cells, Ref. 3^ 
with 2 fiuorite cubic cells, Ref. Ijj 
with 1 fiuorite cubic cell, Ref. 
with 2 fiuorite cubic cells, Ref. l^l 



The PDM equations can be derived by considering only isolated point defect excitations (without internal structure): 
Vb, Vu, lo and lu- In a closed system the particle numbers must be conserved, reduces the number of independent 
defect modes to three. On the other hand, the formation energy reference state for point oxygen and uranium defects 
usually are different, therefore one should instead use three independent combinations of these isolated defects to 
eliminate this ambiguity. The simplest candidates are oxygen and uranium Frenkel pairs and Schottky defect (or 
equivalently, anti-Schottky defect). Consequently, M — 3 and i = 1, • • • ,3 corresponds to the isolated Frenkel pairs 
and Schottky defect, respectively. In this way Eq.® becomes 

PFPo=exp(^^^^^^[Vo][Io], (8) 

PFPu = exp (^^) ^ [Vu][Iu], (9) 



PS = e^p (^-^^ = [VonVu], (10) 

and the composition equation expressed in point defect populations 

x^2{[Vu]-[Iu]) + [Io]-nVo]. (11) 

Eas.(|5K fTT|) comprises the PDM equations exactly. 

To include cluster effects, taking that in C42 configuration as example, we need reinterpret the two intcrstitials 
as an isolated diagonal pair {dp). Assuming this interstitial pair is predominant over the point one, then Eq.([8]) is 
replaced by 

i2 _ ( —Edp — 2Evo 



where two isolated oxygen vacancies have been introduced to eliminate the ambiguity in extrinsic defect formation 
energy. Also, the composition equation becomes 

X = 2 {[Vu] - [lu]) + Pdp - 2[Vo], (13) 

with Eqs.® and PH)) keep unchanged. This procedure can be extended to include other independent clusters easily. 

2. Defect concentrations in PDM 

In the poind defect approximation, the formation energy of a Frenkel pair of the X specie is given by 

Efp,=E^-'+E^'+'-2E^ (14) 
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Composition (-x) Composition (+x) 

FIG. 11: Analysis of the point defect model at a temperature of 1700 K. Variation of the concentrations of point defects with 
the deviation from stoichiometry x: hypostoichiometric regime (on the left) and hyperstoichiometric regime (on the right). 
Solid (respectively dotted and dashed) line indicate the concentration in oxygen interstitial (respectively oxygen vacancy and 
uranium vacancy). The concentration of uranium interstitial is negligible. 



and for the Schottky defect {S) by 

Es^El-' + 2El-'-i^E\ (15) 

with N the number of atoms and the total (or cohesive) energy in the defect-free supercell; Ey^]^ the energy of 

the cell with the defect. Here we use C8±i and "C8±i to model the point defects, thus iV = 96 and E^ and Ey^j^ 
can be obtained by timing 8 to the corresponding cohesive energies listed in table HI 

The formation energies of the defects obtained are listed in table lIVl They are compared to the previous theoretical 
resultsiiii^'i^i^'^ and PDM estimates based on diffusion measurements.^- Note that the GGA+U employed the same 
U parameter as in this work. A detailed comparison of its results with LSDA+U can be found in Ref. z|. Despite it 
produced a similar band gap and local magnetic moment as LSDA+U, it predicted a big lattice constant ~ 5.55 A. By 
figure [7| we know this would lead to an underestimation of the oxygen interaction with the matrix. Table [TVl proves 
this by showing a smaller absolute value of the oxygen interstitial and vacancy formation energies than any other 
calculations. However, this failure is not from GGA but the parameter of U.^^ Besides, this U also underestimates 
the formation energy of uranium interstitial greatly, implying one needs to fitting an own U value for GGA functional 
separately. 

The improvement of LSDA+U over the pure GGA or LDA results is significant. Both the latter underestimate the 
formation energy of uranium vacancy by about 2 times, and 10% 20% for oxygen vacancy. By the lump, LSDA+U 
corrects the energy for 0-FP by a 50% and 38%, a 46% and 61% for U-FP, an 89% and 83% for Schottky defect over 
GGA and LDA, respectively. This correctness makes our LSDA+U results the first ab initio defect formation energies 
that predict the predominance of oxygen defects within a broad enough stoichiometric range over uranium ones (For 
the performances of LDA or GGA formation energies and the PDM anticipation, please see Refs.0Ii3i]). 

The defect concentrations, or equivalently populations, calculated with PDM equations are shown in figure[TTJ An 
arbitrary temperature of 1700 K is chosen. We see oxygen interstitial dominates when x > while it is oxygen vacancy 
when a; < 0. At x ~ 0, 0-FP overwhelms. This picture is in a good agreement with diffusion measurements interpreted 
by PDM,^^ but different from neutron diffraction data where non-negligible oxygen vacancies were observed when 
X > Oi^i^j^jL^^ The population of oxygen vacancy predicted by PDM is too low to be true. To increase this population 
in the regime a; > 0, one needs to take clustering effect into account. 
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Composition (-x) Composition (+x) 

FIG. 12: Analysis of independent clusters model with isolated diagonal oxygen interstitial pair approximation. Others are the 
same as in figure [TT] 

3. Defect concentrations with independent clusters approximation 

Assume the oxygen diagonal pair in CA2 is dominant over the single interstitial, one can calculate the clustering 
effect formally. It is not a promising assumption due to the small energy difference between them. While it can 
be used to analyze the influence of pure interstitial clusters that occupied only the octahedral sites on the vacancy 
populations (they should have similar effects). Also it serves to show how the independent clusters approximation 
works out. 

Using the defect formation energy of CA2 and Eqs.([9l fT0lfT2l and[T3 |) . we calculated the defect populations following 
the same manner as PDM, the result is presented in figure [T^ Note here that ppd turns out to have the same numeric 
value as [Io\- We see that this pure clustering mechanism do decrease the oxygen interstitial population, but that 
of oxygen vacancy in the a; > regime is still too low. Another problem raised here is the population of uranium 
vacancy is closely pinned on that of the oxygen interstitial. It is not what we wanted. Roughly, figure [12] suggests 
clusters associated with oxygen vacancies is necessary in order to enhance the latter 's concentration greatly in a; > 
regime and to pin on that of oxygen interstitial, as implied by the neutron diffraction measurements. 

IV. CONCLUSION 

In summary, we performed a comprehensive calculation on defect properties in U02±x with LSDA+U method. 
The volume changes induced by defects and their formation energies were computed accurately. Analysis of these 
energies for a series configurations concluded that defect clustering is unavoidable when x > 0.03, compatible with 
experimental fact. Atomic charge calculations in Bader's definition, however, showed the difficulty to oxidize uranium 
to and the charged oxygen apt to losing its electrons, against common expectation. As the simplest interstitial 
cluster, oxygen dimer behaviors in a manner similar to a normal oxygen in energetics and charge state. It was identified 
as ionic dioxygen molecule with two excess electrons. Static and vibrational free energy calculations, however, showed 
it quite unstable and might only be a transient state during oxidization process. 

Oxygen dimer is the extreme case for interstitialcy diffusion of oxygen, which may induce a charge fluctuation 
with a magnitude less than 0.2 e. It also presents as a special case for Willis O site occupancy under stretch. The 
stabilization mechanism for this site under ambient conditions, however, is attributed to a V-20 (O ) triple by the 
local stability analysis. Also, O site is stabilized only by the nearest oxygen vacancy pair. This comprises the basic 
clustering pattern for defects in UO2+2:: play with the four building blocks (V-(2)0 and V-(2)0 -V) by sharing the 
vacancy sites. The actual stability of clusters should be judged by the formation energies, which beyond the scope of 
this paper and we prefer to future work. A quasi phase diagram for defect clusters vs composition was also proposed 
to explain the observed population ratios of O and O sites, which of course requires further refinement step by step 
when more calculations and experimental data are available. 
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The formation energy of Frenkel pairs and Schottky defect calculated with LSDA+U have been improved more 
than 50% over the GGA and LDA results. With these energies and the point defect model, we first time showed 
the predominance of oxygen defects by first principles. Finally we generalized the PDM to independent clusters 
approximation that allows us to compute the population of clusters, and revealed the necessity to move on to Willis 
type clusters. 
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